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Lyman-a Emitters in Cosmological Simulations I: Lyman-a 
Escape Fraction and Statistical Properties at z = 3.1 
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ABSTRACT 

We use very large cosmological Smoothed-Particle-Hydrodynamics simulations to 
study the properties of high redshift Lyman-a emitters (LAEs). We identify star- 
forming galaxies at z = 3.1 in a cosmological volume of 100 /i^^Mpc on a side. We 
develop a phenomenological model of absorption, scattering and escape of Lyman-a 
photons on the assumption that the dumpiness of the inter-stellar medium in a galaxy 
is correlate d with the large r scale substructure richness. The radiative transfer effect 
proposed bv lNeufeld] (|1991[ ) allows a large fraction of Lyman-a photons to escape from 
a clumpy galaxy even if it contains a substantial amount of dust. Our model repro- 
duces, for the first time, all of the following observed properties of LAEs at z ~ 3.1: the 
angular correlation function, ultra-violet and Lyman-a luminosity functions, and the 
equivalent width distribution. A simple model that takes only dust absorption into 
account fails in matching the observational data, suggesting that the kind of effect 
we consider is needed. Our model also predicts a bimodal age distribution for LAEs. 
There are old, massive and dusty LAEs, similar to recently found high redshift LAEs. 
The large LAEs have escape fractions of Lyman-a photons of fcsc ^ 0.05 — 0.1. 
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1 INTRODUCTION 

A population of star-forming galaxies at high redshifts 
are characterized by their strong Lyman-a (Lya) line 
emission. Such Lya emitters (LAEs) have been found at 
various redshifts by na rrow-band surveys u sing 8-10 m 



llMori. Umemura. fc Ferrara.1 
l2006l ). 



12004 iMori. fc Umemura.1 



20031: Shimasaku et all 


2003, 


2006: Havashino et al.l 


2004; 


Ouchi et all 2004, 120051. i2008i: iTanieuchi et al.l 


2005; 


Matsuda et al.ll2004l, l2005l: live 


et al.ll2008t). 





It is generally thought that the strong Lya emission 
physically originates from star-forming regions (Hll regions) 
in a young starburst galaxy. Some LAEs have very large 
equivalent widths (EWLya) exceeding 400A, which is diffi- 
cult to explai n with ordinary s tellar p opulation s ynthe sis 
models (e.g., ICharlot. fc Fail] (| 19931 ): ISchaerer.1 ('2003')). 
Alternative physical models include cooling radiation from 
a primordial collapsing gas (|Haiman. Spaans. fc Q.uataert.1 
bOOd: iFardal et aUl200ll). fr om a galactic wind-driven shell 
ijTaniguchi. fc Shiova.ll2000l ), and from supernova remnants 



* E-mail:ikko. shimizu@ipmu.jp 

f E-mail:naoki. yoshida@ipmu.jp 

J E-mail:tokamoto@ccs. tsukuba. ac.jp 



Recent large LAE surveys provided an array of statis- 
tical properties of LAEs such as the Lya luminosity func- 
tion, two-point angular correlation function and the evo- 
lution of them. The observations generally suggest that 
LAEs are not simply a subset of star-forming galaxies. In- 
deed, theoretical models proposed so far do not fully ex- 
plain the observed properties. In particular, reproducing 
very large equivalent widths of some bright LAEs appears 
to be challenging. Lya photons are easily absorbed by dust 
and thus it is naively expected that LAE is a very young 
and dust-free galaxy. While some observati ons and theo- 
retical studies actually support the notion |Gawiser et al.l 
l2006l . l2007l : iMori. fc Umemura.ll2006l : IShimizu et al.ll2007l ). 
more recent multi-wavelength observations of LAEs in op- 
tical, infrared and sub-millimeter suggest that there are 
LAEs tha t are indeed old and dusty ('Finkelstein et al."'200^; 
Lai et all [2008: Matsuda et al, 2007; U chimoto et aL 2008 



Finkelstein et all l2009d : iTamura et all I2OO9I : lOno etld 
2010!) . Interestin gly, such a populat ion increases with de- 
creasing redshift l|Nilsson et al.ll2009l ). There is even an ev- 
idence that some sub-millimeter galaxies show strong Lya 
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emission l|Smail et al.ll200^ ). The existence of a substantial 
amount of dust appears incompatible with strong Lya emis- 
sion. There must be a physical mechanism for Lya photons 
to escape from such dusty galaxies. 

There have been a number of theore t ical s tudies on 
the population of LAEs. iNagamine et al.l (|2010l ) study a 
stochastic m odel where a galaxy goes through LAE phase 
occasionally. IShimizu et all l|2010h argue that an old evolved 
galaxy can become a LAE as a consequence of delayed gas 
accretion. Their model assumes that the delayed starburst 
occurs at outskirt of a galaxy so that Lya photons can 
escape easie r than those emitted from the central region. 
iDaval et al] (|2009l , l20ld . I2OIII ) perform cosmological sim- 
ulations, whi ch reproduce well the UV and Lya luminos- 
ity functions. iDaval et al.l (|201ll ) also calculate the neutral 
hydrogen fraction at 2 = 5.7 by using a combination of 
their LAE formation model and radiation transfer calcu- 
lation of hy d rogen reionization. In the most recent work of 
iDaval et al.l (|201ll ). absorption of Lya photons by dust is 
treated in a very simple manner, where a dusty gas is a 
pure absorber. 

It has been suggested that large-scale gas motions 
in and arou nd galaxies can a ffect the absorption of 
Lya photons (jZheng et al.ll2010al lbh. Strong galactic winds 
from Lya e mitting galaxies are indeed found in the lo- 
cal universe (iLequeux et al.lll995l: iKunth et al.lll998l . l2003l : 



iMas-Hesse et al.l l2003l : iKeel et aLlboOSl ). Strong winds are 
also seen in high-z L yman break galaxie s (LBGs) with 



strong Lya 

'2OO6; Bower et al.l |2004|: 



Pettini et al.l l2002l: [shaplev et al 



Wilman et al 



2005: iFrve et al 



20071 : iPenTericci et al.l l2007l : iTapken et aimOOT ). Although 



the large-scale velocity structure is rather important, it is 
unlikely that all the LAEs blow strong galactic winds be- 
cause only a small fraction of high rcdshift galaxies shows 
the signature of strong outflow (McLindon ct al. 2010). It 
appears that another physical process is necessary for mas- 
sive du^tygalaxies to be bright LAEs. 

iNeufeld] (|l99lh proposed an important effect for Lya 
transfer. In a clumpy, multi-phase inter-stellar medium 
(ISM), dust is locked up in small cold clouds. Lya pho- 
tons can then escape from the clumpy ISM easier than 
continuum photons because Lya photons, having a very 
large scattering cross-section, are preferentially scattered at 
the surface of the clouds. On the other hand, UV contin- 
uum photons are easily absorbed by dust. Neufeld's model 
can explain not only the existence of Lya emission from 
dusty galaxies b u t also t he observed high equivalent widths 
llHansen. fc Oh.l I2OO6I: iKobavashi. Totani fc Nagashima.1 



I2OO7I . I2OIGI : iFinkelstein et al 1 12008| . l2009al lbl). Clearlv it is 
important to incorporate the effect in modelling LAEs. 

In this paper, we study the statistical properties of 
LAEs at z = 3.1. We perform large cosmological hydrody- 
namic simulations for the standard ACDM cosmology. Our 
simulations follow star formation, supernova feedback, and 
metal enrichment. For galaxies identified in our cosmologi- 
cal simulation, we calculate the spectral evolution and dust 
extinction. Unlike in semi-analytic methods, we can directly 
study the internal structure of simulated galaxies as well as 
their spatial distribution in a cosmological volume. 

Throughout this paper, we adopt the ACDM cosmol- 
ogy with the matter density JIm ~ 0.27, the cosmological 
constant f^A ~ 0.73, the Hubble constant /i = 0.7 in units of 



Ho = 100 km s"^ Mpc"\ the baryon density = 0.046, 
and the matter density fluctuation s are normalized by set- 
ting (78 — 0.81 l|Spergel et al.|[2"003l ). All magnitudes are ex- 
pressed in the AB system, and all Lya EWLya values in this 
paper are in the rest frame. 



2 THEORETICAL MODEL 
2.1 Numerical simulations 

Our simulation code is based on an early version of the Tree- 
PM smoothed particle hydrodynamics (SPH) code gadget- 

3 which is a s uccessor of Tree-PM SPH code gadget-2 
(|Springel.ir2005l ). We simulate N = 2 x 640^ particles in 
a comovmg volume of 100 ft"^ Mpc on a side. The mass 
of a dark matter particle and that of a gas particle are 
2.41 X lO*/i"^M0 and 4.95 x lO^ft'^M©, respectively. 

Physical processe s such as star formation and feedbac k 
are implemented as i nlOkamoto. Nemmen fc Bower. I l|2008t ): 
'Oka moto. fc Frenk.l (|2009l ): lOkamoto et al. I 12010'). In par- 
ticular, our simulations employ a new galactic wind model 
in which the initial velocity of a wind particle is proportional 
to the local velocity dispersion of the dark matter parti- 
cles. This is motivated by observations that suggest large 
scale outflows have velocities t hat scale with the circular 
velocity of their host galaxies l|Martin.l |2005| ). As a proxy 
for host halo's circular velocity, which is not easily calcu- 
lated on-the-fly, we use the local one-dimensional velocity 
disper sion, determined from neighbouring dark matter par- 
ticles. lOkamoto et al.l |20l3) found that this quantity, cr, 
is strongly correlated with the maximum circular velocity 
of host (sub-) halos, Wmax, and the relation between these 
quantities does not evolve with redshift. This prescription 
results in a wind speed that increases as a halo grows and 
hence a wind mass-loading (wind mass per unit star forma- 
tion rate) is highest at early times (or in small halos). This 
scaling has been shown to repro duce the physical proper- 
ties of the local grou p satellites (jOkamoto. fc Frenk]|2009l : 
lOkamoto et al.ll2010l ). 

Our simulations includ e the time-evolvi ng pho - 
toionization background (|Haardt. fc Madau.l 1200 ll ). 
metalli city-dependent gas cool ing an d photoheat- 
ing (jWiersma. Schave. fc Smith] |2009| ) , supernovae 
feedb ac k and chemical enrichme n t ([Okamoto et al.l 
I2OO5I : lOkamoto. Nemmen fc Bower.l \2Qm . We use 



metallicity-dependent stellar l i fetim es and ch e mical yields 
llPortinari. Chiosi. fc Bressan.1 Il998l : lMarigo.l [200J). The 
details of these processes are found in the above references. 
Here we give a brief description. Each SPH particle can 
spawn a new star particle when the particle satisfies a 
set of standard criteria for star formation. A star particle 
carries its properties such as mass, formation time and 
metallicity. We calculate the spectral energy distribution 
(SED) of each star particle using the populatio n synthesis 
code pegase l|Fioc fc Rocca-Volmerange.l 1 19971 ). We then 
sum the individual SEDs to obtain the total SED of a 

simulated galaxy. 

W e run a friends-of-friends group finder (iDavis et al.l 
Il985h to locate groups of stars, i.e., galaxies. We also iden- 
tify substructures (subhalos) in each FoF group us ing SUB- 
FIND algorithm developed bv lSpringel et al.l (|200ll ). For the 
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identified galaxies, we calculate the intrinsic Lya luminosi- 
ties using 'PEGASE'. In PEGASE, two thirds of ionizing 
photons are converted to Lya photons under the assump- 
tion of the case B recombination. Finally, we calculate the 
effect of dust extinction on ultra-violet continuum and Lya 
emission. This is the key component of our model. We will 
describe the details in the next subsection. 



2.2 Model description for dust absorption 

We calculate the SED and the Lya luminosity for individual 
galaxies identified in our cosmological simulation. In order to 
compare our model predictions directly with observational 
data, we need to include the effect of dust absorption. We 
first assume that the dust mass to metal mass ratio is a 
constant of 0.4, consistent with the local value. Note that 
the metallicity of a galaxy is calculated from the metallic- 
ities of gas particles. We further assume that the galaxies 
are roughly spherical; we use only one length scale, the ef- 
fective radius, to evaluate the optical depth. We calculate 
the optical depth rd(A) for UV continuum photons as 

3Ed 



(1) 



where ad and s are the typical size of dust grains, and 
the material density of dust grains, respectively. We adopt 
the standard choice of Od = 0.1 and s = 2.5 g cm"'^ 
iTodini fc Ferarra.ll200il : lNozawa et al.ll2003l '). The dust sur- 
face mass density Ed is 

Md 



Erf = 



(2) 



where Md and rd are the total dust mass (40% of metal 
mass) and the effective radius of the galaxy, respectively. 
The effective radius rd is a fraction of the virial radius and 
is given by /d^vir with /d = 0.18. The escape fraction of UV 
continuum photons /cont(A) is then calculated as 



/cont(A) — 



exp(-rd(A)) 
^d(A) 



(3) 



We use the du st optical constant Q as a function of wave- 
length given inl Praine fc Lee" We also calculate the 
IGM a bsorption at the blue side of 1216A following iMadau"] 
ll 199511 . 

A crucial quantity in our model is the effective optical 
depth of Lya line, TLya . We employ two models to calculate 
the effective optical depth. A simplest assumption would be 
to set 



abs 



(4) 



where Cabs is a constant parameter. Then, the escape fraction 
of Lya photons fi^y% is given by 



/abs 
Lya 



exp ( 



abs \ 



(5) 



'Lya 



We call this model as pure absorption model. Essentially, 
dust is treat ed as a pur e absorber of Lya photons in this 
model (e.g.. iDaval et al.1 ([2009, 2010. 1201 il l). Lya photons 
can not easily escape from dusty galaxies which have a large 
effective optical depth. 

The second model, which we propose in the present pa- 
per, is motivated by the multiphase inter-stellar medium 



(ISM) model of lNeufeld.1 l|l99lh . In a clumpy ISM, Lya pho- 
tons are scattered mostly at the surface of cold clumps before 
they are absorbed by dust. A large fraction of Lya photons 
can then escape from a clumpy ISM through multiple scat- 
ters. The key quantity here is the overall dumpiness of the 
ISM in a galaxy. Since our cosmological simulation does not 
resolve the fine structure of the ISM, we need to estimate the 
dumpiness of the ISM in some way. It is probably reason- 
able to expect that the ISM structure is well developed in a 
massive galaxy which itself has rich substructures. We make 
an assumption that the dumpiness of the ISM has a cor- 
relation with the larger scale internal structure of a galaxy. 
In practice, we use the number of subhalos (satellites) A'^sub 
of the galaxy as a measure of its "dumpiness". We intro- 
duce the dumpiness factor S, which depends on A'^ub, and 
express the effective optical depth of Lya line as 



sub Q 

'^LyQ Cgub'^'^di 



(6) 



where Csub is a normalization constant. For simplicity, we 
assume the dumpiness factor is expressed as a power-law 



(7) 



After some experiments, we find that setting a = —0.5 works 
remarkably well, reproducing nearly all the important obser- 
vational data, as will be shown in the following sections. We 
call this model the substructure model. The escape fraction 
of Lya photons /^ya in this model is given by 



/■sub 

J Lya — 



l-exp(-r£fj 



(8) 



'Lya 



The combination of equation ((7)l and equation (|8]) ef- 
fectively yields a larger escape fraction for a dusty but 
"clumpy" galaxy. Finally, we fix the normalization constants 
in the two models, Cabs and Csub, by matching the num- 
ber density of simulated LAEs wit h the observed nu mber 
density jilaes 5 x 10"^ [Mpc'^] l|Ouchi et al.ll2008D . We 
identify simulated galaxies which satisfy EWLya > 20A and 
Ll];^ > 1.0 X W^^lerg/s] as LAEs. 



3 RESULT 

3.1 The Spatial Distribution of Simulated LAEs 

Fig. [T] shows the projected distribution of the LAEs at z = 
3.1 in the simulation volume of 100 comoving h~^Mpc on a 
side. The point size represents the Lya luminosity of each 
galaxy. Luminous LAEs are shown by large and red points. 
The overall distribution appears quite similar in the two 
models, with LAEs approximately tracing the underlying 
matter distribution. The luminous LAEs are clustered more 
strongly in the substructure model, although the difference 
is small. 

We quantify the spatial distribution of the simulated 
LAEs. We calculate the two-point angular correlation func- 
tion (ACF) and compare it with the observed correlation 
function. The AGE o f LAEs found in the SSA22 field 
(|Havashino et al1l2004l ) is used for the comparison. Fig. [2] 
shows the ACFs for the two models. Both of the models re- 
produce the observational result well. The LAEs are hosted 
by dark halos with mass of ~ 1O"M0. 
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The pure absorption model 



The substructure model 
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Figure 1. The spatial distribution of simulated LAEs. Tiie left and right panels show LAEs in the pure absorption model and those in 
the substructure model, respectively. The point size is scaled with Lya luminosity of each galaxy so that luminous galaxies appear as 
large points. The smallest and biggest points correspond to LAEs with lO**^ ^Lya lO'^^'^ergs"^ and I/Lya ^ 10'*''ergs~^, respectively. 
The red-points shows brightest LAEs (> 6 X 10*^ [erg/s]). 




10 20 
correlation length [Mpc] 

Figure 2. The two-point angular correlation function (ACF). 
The solid and dashed lines represent the substructure model and 
the pure absorption model, respectively. Point s with error bar 
are t he ACF of LAEs observed in SSA22a field l|Havashino et al.l 
|2004) . 



3.2 Lya Luminosity Functions 

In Fig. |3l we compare the Lya lu minosity functions witii 
tiie observational data at 2 = 3.1 (jOuclii et al.ll2008l ). Tlie 
pure absorption model (left panel) does not reproduce tiie 
observational data at tiie briglit end. It predicts a smaller 
number of LAEs by a factor of ten at LLya ~ lO^^ergs"^. 
In botli our models, tiie intrinsic Lya luminosity of a galaxy 
is proportional to its star formation rate (SFR), and tlius 



tlie intrinsically brigiit LAEs are liosted typically by mas- 
sive lialos. However, such massive galaxies are also aged and 
dusty. Because of the strong dust absorption of Lya photons 
(see equation [4]), dusty star- forming galaxies do not appear 
as LAEs in the pure absorption model. 

Our substructure model, on the other hand, shows a 
substantially better agreement with the observational data. 
In the substructure model, dusty, aged and massive galax- 
ies have complex internal structures. We assume that the 
evolved galaxies have also a more complex ISM structure, 
where Lya photons can escape easily because of the Neufeld 
effect. Thus, even aged and dusty galaxies appear as bright 
LAEs and the resulting luminosity function matches the ob- 
servation very well at the bright-end. 

Interestingly, the substructure model predicts that the 
Lya escape fraction of massive galaxies is /oac ~ 0.05 — 0.1, 
with a substantial dispersion. In the pure absorption model, 
the Lya escape fraction decreases proportionally to the halo 
mass because of the size effect and also because massive 
galaxies are more metal-enriched. 



3.3 Mt 



EWt 



distribution 



Fig. U shows the distributions of our simulated LAEs in the 
Muv — EWLya plane ai z = 3.1. In this plot, both models 
appear similar, with the LAEs populating the region where 
the observational data points also lie. Although both the 
models reproduce the observed trend that the equivalent 
width decreases with UV luminosity, only the substructure 
model predicts the existence of UV bright galaxies (Muv < 
—20) with large EWLya > lOOA. This is again owing to the 
effect of enhanced Lya luminosity in the substructure model. 
We also note here that the substantial dispersion of EWLya 
for UV bright LAEs, as is also found in the observations, is 
important. Apparently the Lymana escape fraction is not 
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Figure 3. The Lyo luminosity function for the pure absorpt i on mo del (left) and for the substructure model (right). We compare the 
model predictions with the observational data of lOuchi et al.l (|2008l ) with error bars. The filled circle marked an open square indicates 



LAE with AGN. 





Figure 4. Distributions of LAEs i n the M^jy — EWLv g plane at z = 3.1. The simulated LAEs are shown by dots, whereas points with 
error bars are observational data of lOuchi et al ] ll2008r ). The filled circles marked an open square indicate LAEs with AGN. 

a simple function of the galaxy mass only. We discuss this tinuum is absorbed by dust significantly but Lya photons 
issue more in detail in the discussion section. can escape again by the scattering effect. 



3.4 Lya Equivalent Width Distribution 

In Fig. O we compare the rest-frame Lya equivalent width 
of simulated LAEs (histograms) with the observational 
data (solid line with error bars). In the pure absorption 
model, EWLya distribution is approximately-constant at 
EWlyq > lOOA whereas the observed distribution decreases 
toward high EWhya ■ The substructure model reproduces the 
observed EWLya distribution very well. It is interesting that 
the fraction of large EWLya galaxies is similar between the 
two models. However, the substructure model predicts more 
LAEs with EWLyc ~ 100 - 200A. For these LAEs, UV con- 



3.5 Age distribution of simulated LAEs 

Finally, we show the age distribution of simulated LAEs in 
Fig. [6l We assign the age of a simulated LAE by calculat- 
ing the mass weighted mean of the stellar ages. The simu- 
lated LAEs shows clearly a bimodal distribution. Recently 
formed galaxies are young LAEs, whereas massive, dusty 
galaxies with tage > 10* years can also appear as bright 
LAEs. The latte r popul ation is indeed found recently by 
iFinkelstein et all l|2009bl ). 
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Figure 6. The age distribution of simulated LAEs for the sub- 
structure modeh There are two types LAEs; young LAEs that are 
populated at the left of the plot and old LAEs that arc broadly 



distributed around ta 



' 2 X 10° years. 



4 CONCLUSIONS AND DISCUSSION 

We study a number of observed properties of LAEs us- 
ing a large cosmological hydrodynamic simulation. The 
simulation follows the formation and evolution of star- 
forming galaxies by employi n g new feedback models of 
lOkamoto. Nemmen fc Bower .1 ( 20081 ) . We develop a novel 
model in which Lya photons can escape from a massive, 
dusty galaxy if the ISM is clumpy. The idea is motivated by 
the multiple scatterin g and e scape of Lya photons originally 
proposed bv lNeufeld.1 (|l99lh . 

Our physical model of LAEs reproduces not only the 
Lya luminosity function but also Muv — EWlyq distribu- 



tion, the equivalent width distribution and the angular two- 
point correlation function at z — 3.1. It is the first theo- 
retical model that successfully reproduces all these observa- 
tional results using cosmological hydrodynamic simulations. 
Interestingly, the model predicts the Lya escape fraction is 
roughly constant for the massive LAEs. Massive, dusty and 
aged galaxies can appear as bright LAEs, as is consistent 
with recent observations. Contrastingly, the pure absorption 
model fails in reproducing many of the observational results. 
We argue that the kind of effect to enhance Lya photon es- 
cape from massive dusty galaxies is necessary for modeling 
LAEs. 

It is interesting to ask if the Lya escape fraction can 
be described as a simpler function of galaxy mass. In our 
simulation, the substructure abundance of simulated galax- 
ies roughly scales with their host halo mass, although with 
substantial dispersions. It has turned out that the dispersion 
is also important to reproduce some observational results. To 
explicitly test the idea, we have employed a simpler model 
in which the Lya escape fraction is a function of host halo 
mass so that the escape fraction approximately matches to 
the mean escape fraction of our substructure model for a 
given halo mass. We have found that this simple model fails 
in reproducing the Muv — EWLya distribution and ACF. 
The simple model boosts the Lya luminosity too much and 
cannot explain the existence of UV-bright LAEs with low 
equivalent widths (see Fig. 3]). Furthermore, in the simple 
model, massive galaxies are identified selectively as LAEs. 
Consequently the strength of ACF becomes large, although 
it is still within the upper bound of the current data. Over- 
all, we conclude that the internal structure of a galaxy is 
an important factor to determine it s appearance as a LAE. 
Interestingly, ICharlot" fc Fall.1 (1199311 argue that a large scat- 
ter of Lya emission versus metallicity correlation might be 
caused by the structure of the interstellar medium. 

In the present paper, we have focused on the physi- 
cal properties of LAEs at a particular redshift of z = 3.1 
where an array of observations of LAEs are available. Ac- 
cording to very recent observations, the fraction of old pop- 



Simulating Lya emitters 7 



ulation in observed LAEs increases with decreasing redshift, 
whereas the escap e fraction of Lyg emission becomes lar ger 
at higher redshift (jNilsson et al.ll2009l : iHaves et al.ll2010l l. It 
is certainly interesting and important to test whether or not 
our model can reproduce not only the physical properties 
of observed LAEs but also their evolution. We will present 
the properties of simulated LAEs at various redshifts in a 
forthcoming paper (Shimizu et al. 2011 in prep). We will 
also study the epoch of hydrogen reionization using our the- 
oretical model. 
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